{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "#from vk_download import main\n",
    "import scipy.sparse as sprs\n",
    "import vk\n",
    "import time\n",
    "from tqdm.auto import tqdm\n",
    "import os\n",
    "import gc\n",
    "import numpy as np\n",
    "from scipy.sparse import csr_matrix, lil_matrix\n",
    "import matplotlib.pyplot as plt\n",
    "import networkx as nx\n",
    "import seaborn as sns\n",
    "import scipy\n",
    "import random\n",
    "from scipy import stats\n",
    "import chart_studio.plotly as py\n",
    "import plotly.graph_objs as go\n",
    "from scipy.stats import binom\n",
    "import math"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#SUPPRESS SCIENTIFIC NOTATION\n",
    "pd.set_option('display.float_format', lambda x: '%.3f' % x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sub_hist(data):\n",
    "    bins = np.logspace(0, 3.4, 20)\n",
    "    widths = (bins[1:] - bins[:-1])\n",
    "    hist = np.histogram(data, bins=bins)\n",
    "    hist_norm = hist[0]/widths\n",
    "    plt.figure(figsize=(5, 5))\n",
    "    plt.rcParams.update({'font.size': 15})\n",
    "    plt.bar(bins[:-1], hist_norm, widths, color='grey')\n",
    "    plt.xscale('log')\n",
    "    plt.yscale('log')\n",
    "    plt.xlabel('Число друзей $k$', fontsize = 10)\n",
    "    plt.ylabel('Число пользователей, имеющих $k$ друзей', fontsize = 10)\n",
    "    plt.show()\n",
    "    return 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Download data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1660927\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>$p_{i}(t_{1})$</th>\n",
       "      <th>$p_{i}(t_{2})$</th>\n",
       "      <th>$p_{i}(t_{3})$</th>\n",
       "      <th>fr_number</th>\n",
       "      <th>$p_{-i}(t_{1})$</th>\n",
       "      <th>$p_{-i}(t_{2})$</th>\n",
       "      <th>$p_{-i}(t_{3})$</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.515</td>\n",
       "      <td>0.515</td>\n",
       "      <td>0.512</td>\n",
       "      <td>7.000</td>\n",
       "      <td>0.410</td>\n",
       "      <td>0.467</td>\n",
       "      <td>0.457</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.186</td>\n",
       "      <td>0.186</td>\n",
       "      <td>0.186</td>\n",
       "      <td>7.000</td>\n",
       "      <td>0.545</td>\n",
       "      <td>0.577</td>\n",
       "      <td>0.578</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.604</td>\n",
       "      <td>0.603</td>\n",
       "      <td>0.604</td>\n",
       "      <td>6.000</td>\n",
       "      <td>0.334</td>\n",
       "      <td>0.309</td>\n",
       "      <td>0.287</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.608</td>\n",
       "      <td>0.615</td>\n",
       "      <td>0.638</td>\n",
       "      <td>3.000</td>\n",
       "      <td>0.411</td>\n",
       "      <td>0.414</td>\n",
       "      <td>0.413</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0.730</td>\n",
       "      <td>0.723</td>\n",
       "      <td>0.727</td>\n",
       "      <td>3.000</td>\n",
       "      <td>0.605</td>\n",
       "      <td>0.616</td>\n",
       "      <td>0.615</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   $p_{i}(t_{1})$  $p_{i}(t_{2})$  $p_{i}(t_{3})$  fr_number  $p_{-i}(t_{1})$  \\\n",
       "0           0.515           0.515           0.512      7.000            0.410   \n",
       "1           0.186           0.186           0.186      7.000            0.545   \n",
       "2           0.604           0.603           0.604      6.000            0.334   \n",
       "3           0.608           0.615           0.638      3.000            0.411   \n",
       "4           0.730           0.723           0.727      3.000            0.605   \n",
       "\n",
       "   $p_{-i}(t_{2})$  $p_{-i}(t_{3})$  \n",
       "0            0.467            0.457  \n",
       "1            0.577            0.578  \n",
       "2            0.309            0.287  \n",
       "3            0.414            0.413  \n",
       "4            0.616            0.615  "
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data = pd.read_csv('X_opinions.csv')\n",
    "\n",
    "data.drop(['Unnamed: 0'], axis=1, inplace=True)\n",
    "\n",
    "print(data.shape[0])\n",
    "\n",
    "data.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Users' opinions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.48 0.177\n",
      "0.477 0.179\n",
      "0.474 0.18\n"
     ]
    }
   ],
   "source": [
    "print(round(data['$p_{i}(t_{1})$'].mean(), 3), round(data['$p_{i}(t_{1})$'].std(), 3))\n",
    "print(round(data['$p_{i}(t_{2})$'].mean(), 3), round(data['$p_{i}(t_{2})$'].std(), 3))\n",
    "print(round(data['$p_{i}(t_{3})$'].mean(), 3), round(data['$p_{i}(t_{3})$'].std(), 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, '$p_{i}(t_{3})$')"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAFgCAYAAACmDI9oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAXA0lEQVR4nO3df6xkZ3kf8O8T20AbCKZ4VSxjWMqPkEC8/FgZKBIioEg2jeyqQGsUBUxBWxGoSUWjkvwBLVVVJWlJRJxiLGMZKAUiE6ENcoIQUAxNIawd22CMqxX9wbpuvNjBhkBMlz79447h9eWu7+zuzJ1z7/18pNGdOeedM491/Ox8550z51R3BwAAWPMTqy4AAACmREAGAICBgAwAAAMBGQAABgIyAAAMBGQAABisNCBX1dVVdVdVfWWOsZdW1dGquml2e/1W1AgAwO6y6hnka5JccALjP9Ldz5rdrlpSTQAA7GIrDcjdfX2Se8ZlVfXkqvqTqrqhqj5XVU9fUXkAAOxCq55B3siVSf5pdz83yT9P8h+GdS+vqluq6tqqOnc15QEAsJPVqi81XVV7k3y8u59ZVY9McjTJ7cOQh3f3z1TVY5N8p7vvr6p/kuQfdfdLtr5iAAB2sqkF5J9Kcnt3n73Jc05Lck93P3oLSgQAYBeZ1CEW3X1fkv9eVa9Mklqzb3Z/DM0XJbltBSUCALDDrfo0bx9K8l+T/HRVHamq1yX5pSSvq6qbk9ya5OLZ8Muq6tbZ8suSXLqKmgEA2NlWfogFAABMyaQOsQAAgFU7fVUvfNZZZ/XevXtX9fKwo91www3f7O49J/o8fQnLoy9heo7XlysLyHv37s2hQ4dW9fKwo1XV/zyZ5+lLWB59CdNzvL50iAUAAAwEZAAAGAjIAAAwEJABAGAgIAMAwEBABgCAgYAMAAADARkAAAYCMgAADARkAAAYCMgAADAQkAEAYHD6qgtgtapqrnHdveRKAACmwQwyAAAMBGQAABgIyAAAMBCQAQBg4Ed6ABPjx7MAq2UGGQAABgIyAAAMBGQAABgIyAAAMPAjPQCATfjx7O5iBhkAAAYCMgAADARkAAAYbBqQq+oRVfVnVXVzVd1aVf9qgzEPr6qPVNXhqvpiVe1dRrEAALBs88wg35/kJd29L8mzklxQVc9fN+Z1Sf6yu5+S5HeS/OZiywQAgK2xaUDuNd+ZPTxjdlv/E82Lk7xvdv/aJC+teX/uCQAAEzLXMchVdVpV3ZTkriSf7O4vrhtyTpJvJEl3H0tyb5LHbrCdA1V1qKoOHT169NQqBxZCX8L06EtYrbkCcnf/oLufleTxSc6vqmeezIt195Xdvb+79+/Zs+dkNgEsmL6E6dGXsFondBaL7v5Wks8kuWDdqjuSnJskVXV6kkcnuXsRBQIAwFaa5ywWe6rqzNn9v5HkF5J8bd2wg0leM7v/iiSfbpeSAQBgG5rnUtNnJ3lfVZ2WtUD9B9398ap6R5JD3X0wyXuTfKCqDie5J8klS6sYAACWaNOA3N23JHn2BsvfNtz/6ySvXGxpAACw9VxJDwAABgIyAAAMBGQAABgIyAAAMBCQAQBgICADAMBAQAYAgIGADAAAAwEZAAAGAjIAAAwEZAAAGAjIAAAwEJABAGAgIAMAwEBABgCAgYAMAAADARkAAAYCMgAADARkAAAYCMgAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMBCQAQBgICADAMBAQAYAgIGADAAAAwEZAAAGAjIAAAw2DchVdW5VfaaqvlpVt1bVmzcY8+Kqureqbprd3raccgEAYLlOn2PMsSRv6e4bq+pRSW6oqk9291fXjftcd//i4ksEAICts+kMcnff2d03zu5/O8ltSc5ZdmEAALAKJ3QMclXtTfLsJF/cYPULqurmqvrjqnrGcZ5/oKoOVdWho0ePnnCxwOLpS5gefQmrNXdArqpHJvlokl/t7vvWrb4xyRO7e1+S30vysY220d1Xdvf+7t6/Z8+ek60ZWCB9CdOjL2G15grIVXVG1sLxB7v7D9ev7+77uvs7s/vXJTmjqs5aaKUAALAF5jmLRSV5b5LbuvudxxnzuNm4VNX5s+3evchCAQBgK8xzFosXJvnlJF+uqptmy34jyROSpLuvSPKKJG+oqmNJvpfkku7uJdQLAABLtWlA7u7PJ6lNxlye5PJFFQUAAKviSnoAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMBCQAQBgICADAMDg9FUXwPZQVXON6+4lVwIAsFxmkAEAYCAgAwDAQEAGAICBgAwAAAM/0gPYpvx4FmA5zCADAMBAQAYAgIGADAAAA8cgAwAsiN8G7AxmkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMBCQAQBgICADAMBAQAYAgIGADAAAg00DclWdW1WfqaqvVtWtVfXmDcZUVb2rqg5X1S1V9ZzllAsAAMt1+hxjjiV5S3ffWFWPSnJDVX2yu786jLkwyVNnt+cleffsLwAAbCubziB3953dfePs/reT3JbknHXDLk7y/l7zhSRnVtXZC68WAGCBqmquG7vLCR2DXFV7kzw7yRfXrTonyTeGx0fy4yEaAAAmb+6AXFWPTPLRJL/a3fedzItV1YGqOlRVh44ePXoymwAWTF9uHTNVzEtfwmrNFZCr6oysheMPdvcfbjDkjiTnDo8fP1v2IN19ZXfv7+79e/bsOZl6gQXTlzA9+hJWa56zWFSS9ya5rbvfeZxhB5O8enY2i+cnube771xgnQAAsCXmOYvFC5P8cpIvV9VNs2W/keQJSdLdVyS5LsnLkhxO8t0kr118qZwIX9MCAJycTQNyd38+yUOmre7uJG9cVFEAALAqrqQHAAADARkAAAYCMgAADARkAAAYCMgAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMBCQAQBgICADAMBAQAYAgIGADAAAAwEZAAAGAjIAAAwEZAAAGAjIAAAwEJABAGAgIAMAwEBABgCAgYAMAAADARkAAAYCMgAADARkAAAYCMgAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABpsG5Kq6uqruqqqvHGf9i6vq3qq6aXZ72+LLBACArXH6HGOuSXJ5kvc/xJjPdfcvLqQiAABYoU1nkLv7+iT3bEEtAACwcos6BvkFVXVzVf1xVT3jeIOq6kBVHaqqQ0ePHl3QSwOnQl/C9OhLWK1FBOQbkzyxu/cl+b0kHzvewO6+srv3d/f+PXv2LOClgVOlL2F69CWs1ikH5O6+r7u/M7t/XZIzquqsU64MAABW4JQDclU9rqpqdv/82TbvPtXtAgDAKmx6Fouq+lCSFyc5q6qOJHl7kjOSpLuvSPKKJG+oqmNJvpfkku7upVUMAABLtGlA7u5XbbL+8qydBg4y+zJhUz5DAQBT5Up6AAAwmOdCIQBsY/N8s+NbHYAfMYMMAAADARkAAAYOsQAA2GIOfZo2M8gAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMHAlPQBgx5nnSnVwPGaQAQBgICADAMDAIRYAp8hXuQA7ixlkAAAYCMgAADBwiMU25OtcAIDlMYMMAAADARkAAAYCMgAADARkAAAYCMgAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABgIyAAAMBGQAABgIyAAAMNg0IFfV1VV1V1V95Tjrq6reVVWHq+qWqnrO4ssEAICtMc8M8jVJLniI9RcmeersdiDJu0+9LAAAWI1NA3J3X5/knocYcnGS9/eaLyQ5s6rOXlSB7ExVtekNAGAVFnEM8jlJvjE8PjJb9mOq6kBVHaqqQ0ePHl3ASwOnSl/C9OhLWK0t/ZFed1/Z3fu7e/+ePXu28qWB49CXJPN9q+Obna2jL2G1Tl/ANu5Icu7w+PGzZQAAnKR5P5R295Ir2X0WMYN8MMmrZ2ezeH6Se7v7zgVsFwAAttymM8hV9aEkL05yVlUdSfL2JGckSXdfkeS6JC9LcjjJd5O8dlnFAgDAsm0akLv7VZus7yRvXFhFAACwQq6kBwAAAwEZAAAGAjIAAAwEZAAAGAjIAAAwWMSFQgAAtoyrOrJsZpABAGBgBhngIZipAth9zCADAMDADPKEmKkCAFg9M8gAADAQkAEAYCAgAwDAQEAGAICBgAwAAAMBGQAABk7zxmTNe9q77l5yJQDAbiIgAzA3H1yB3UBABgDYxnxwXTzHIAMAwEBABgCAgYAMAAADARkAAAYCMgAADARkAAAYCMgAADBwHmQAYBLmPZ8vLJsZZAAAGJhBBnYlM1UAHI8ZZAAAGJhB3gJmqgAAtg8zyAAAMDCDzLY37wx9dy+5EuAB+hLYzgRkAIBdwAfX+c11iEVVXVBVt1fV4ap66wbrL62qo1V10+z2+sWXCgAAy7fpDHJVnZbk95P8QpIjSb5UVQe7+6vrhn6ku9+0hBoBAGDLzDODfH6Sw9399e7+fpIPJ7l4uWUBAMBqzBOQz0nyjeHxkdmy9V5eVbdU1bVVde5GG6qqA1V1qKoOHT169CTKBRZNX8L06EtYrUWd5u2Pkuzt7vOSfDLJ+zYa1N1Xdvf+7t6/Z8+eBb00cCr0JUyPvoTVmicg35FknBF+/GzZD3X33d19/+zhVUmeu5jyAABga80TkL+U5KlV9aSqeliSS5IcHAdU1dnDw4uS3La4EgGA7ayq5rrBVGx6FovuPlZVb0ryiSSnJbm6u2+tqnckOdTdB5NcVlUXJTmW5J4kly6xZgAAWJq5LhTS3dcluW7dsrcN9389ya8vtjSAE2cWCoBT5Up6p8Ab8fbiCkIwPfoSmCIBGQCAH/LBdXGneQMAgB3BDDIA7ABm/WBxBGTghHkjBmAnm3xA9kYMAMBWmnxABgCmydmc2KkEZGBb8EYMwFYRkDfgjRgAYPcSkGEdx70DwO4mIAMwefN8cPWhFbbWTu5LARkAeBCHGrLbCcjASnkjBmBqXGoaAAAGu24G2WwVAAAPZdcFZADYzUwUweYEZGBpvBEDsB3tmIDsjRgAgEXYMQEZttpOPv8jbEcu8gPTs1370lksAABgICADAMBAQAYAgIGADAAAAwEZAAAGAjIAAAyc5g2WaLue3gZ2Mn0J0zO1vjSDDAAAAwEZAAAGAjIAAAwEZAAAGAjIAAAwEJABAGAgIAMAwMB5kGECpnb+R0BfwhRtVV/ONYNcVRdU1e1Vdbiq3rrB+odX1Udm679YVXtPqSoAAFiRTQNyVZ2W5PeTXJjkZ5O8qqp+dt2w1yX5y+5+SpLfSfKbiy4UAAC2wjwzyOcnOdzdX+/u7yf5cJKL1425OMn7ZvevTfLSmncOHAAAJmSeY5DPSfKN4fGRJM873pjuPlZV9yZ5bJJvjoOq6kCSA7OH36mq2+d4/bPWb4eVsB8moKrm3Q9PPIFt6svty36YAH3JOvbDBJxqX27pj/S6+8okV57Ic6rqUHfvX1JJzMl+mIZl7Ad9uX3ZD9OgLxnZD9NwqvthnkMs7khy7vD48bNlG46pqtOTPDrJ3SdbFAAArMo8AflLSZ5aVU+qqocluSTJwXVjDiZ5zez+K5J8up33BgCAbWjTQyxmxxS/KcknkpyW5OruvrWq3pHkUHcfTPLeJB+oqsNJ7slaiF6UE/qKiaWxH6ZhKvthKnXsdvbDNExlP0yljt3OfpiGU9oPZaIXAAB+xKWmAQBgICADAMBgMgHZ5aynYY79cGlVHa2qm2a316+izp2sqq6uqruq6ivHWV9V9a7ZPrqlqp6zxFr05QToy9XTl6ynL1dvmX05iYDsctbTMOd+SJKPdPezZrertrTI3eGaJBc8xPoLkzx1djuQ5N3LKEJfToO+nIxroi+Z0ZeTcU2W1JeTCMhxOeupmGc/sGTdfX3WzgZzPBcneX+v+UKSM6vq7CWUoi+nQV9OgL5kHX05Acvsy6kE5I0uZ33O8cZ097EkD1zOmsWZZz8kyctnX1VcW1XnbrCe5Zp3P23F6+jL5dOX24O+3F305fZw0n05lYDM9vFHSfZ293lJPpkfzVIAq6MvYXr05TY2lYDsctbTsOl+6O67u/v+2cOrkjx3i2rjR+bpl616HX25fPpye9CXu4u+3B5Oui+nEpBdznoaNt0P647duSjJbVtYH2sOJnn17Ne5z09yb3ffuYTX0ZfToC+3B325u+jL7eGk+3LTS01vhQlczprMvR8uq6qLkhzL2n64dGUF71BV9aEkL05yVlUdSfL2JGckSXdfkeS6JC9LcjjJd5O8dhl16Mtp0JfToC8Z6ctpWGZfutQ0AAAMpnKIBQAATIKADAAAAwEZAAAGAjIAAAwEZAAAGAjIAAAwEJABAGAgILOh2VVnbqqq16xb/rSq+pezKwfNs53Lq+q9y6kSdhd9CdOjL3cmAZnj+YdJ/laS/7Ru+YVJfqW7vz/ndv5dkl+qqqcssjjYpfQlTI++3IEEZI7nsiQf6O7/u275eUm+Mu9Guvt/JPl8kjcsrjTYtfQlTI++3IEE5F2iqu6pqn9dVR+sqrtnj99TVQ/fYOxTkvzdJNeuW35Xkn+c5Oerqme3l8zx8h/N2qdi/7/BQF/C9OhLEgF5V6iqJyR5TJJ/luSeJK9M8rtJDiR54wZPeWmSv0py87CNSnJRkr9O8m+SvGB2+9M5SvjTJH87yc+d9H8E7DD6EqZHX/KA01ddAFti3+zv5d391tn9T1fVi5JckOSd68Y/N8lt3f3/HljQ3V1V/yfJI5J8vLu/cAKvf2uSHyQ5P8M/IrDL6UuYHn1JEjPIu8W+JN9N8tvrlv+3rP2wIFX19Kq6arb8cUm+ucF2zkvSSb48Lqyqd1fVHVXVG714dx9L8q3ZdoE1S+vLqjq3qj5VVbdV1a1V9VuzWa0f0pewoWW/X362qm6uqluq6tqq+qlxvb6cDgF5d9iX5Pruvnvd8scl+d9J0t1f6+7Xz5Y/Isn9G2znvCRf7+6/Wrf8Q0mes0kN98+2C6xZZl8eS/Ivuvtnkjw7yfOS/IMNnqsv4cGW/X55UXfv6+7zkvyvJL+2wXP15QQIyLvDvswa+wFV9agkL0nyqdnjK6rqktnqe5KcucF2zssGX/l09/Xd/Reb1HDmbLvAmqX1ZXff2d2HZve/n+SWJOdu8Fx9CQ+27PfLe2fb+IkkP5m1Web19OUECMg7XFX9ZJInz26jy2Z/3zf7+6wkN83u357kSRts7ulJvnYSNexJ8jez9hUV7Hpb2ZdV9dgkfz/JJ9Yt15cw2Kq+rKrrkvxFkp9O8lvr1unLiRCQd76fy9p+fnJVva2qfr6q/m2Styc50N3fmn2S/Tv5UUP+lyRPmDXq6L4kL6qqF1XV89cf0/gQ9mftU/I8v+CF3WBL+nJ2Wqprk/xud9+27nn6Eh5sS/qyu1+WtUM2/izJr6x7nr6cCAF559uXtUb9e0kuTvIns/uv7O4/mI15WtaOlXrgV7j/OWtf71ywbltvzdqPFD6V5GPdveGP8jZwQZLPbnBMF+xWS+/LqjotyQeT/Hl3//sNatCX8GBb9n7Z3T/I2oz0q9c9T19OhNO87Xz7kny5u2/J2uloNjJ+XZTu/n5V/ccklyT5wLD880mecSIvPnuTfnnW/rEA1mxFX74nybeTvGX9Cn0JG1pqX1bVY5I8bPjNzsszXGlPX06LGeSdb182P5figxp+5rezdgWgp232AlV1VVUdmd0/Mpz+Jlk7yfr3knx4/pJhx1tqX1bVC5O8Lmtf1/55Vd1UVZcNQ/Ql/Lhlv18+JsnHZ6d4+3KSZyZ587BeX05Izf8tOdvN7Jine5P8Wne/5ySef0mSO7v7s6dQw6uS3NHd15/sNmAn0ZcwPfqS9QRkAAAYOMQCAAAGAjIAAAwEZAAAGAjIAAAwEJABAGAgIAMAwEBABgCAwf8H/SC67AmHAIUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x360 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axs = plt.subplots(1, 3, sharey=True, tight_layout=True, figsize=(10, 5))\n",
    "labelsize = 15\n",
    "plt.rcParams['xtick.labelsize'] = labelsize\n",
    "plt.rcParams['ytick.labelsize'] = labelsize\n",
    "plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))\n",
    "n_bins = 20\n",
    "axs[0].hist(data['$p_{i}(t_{1})$'], color='k', bins=n_bins)\n",
    "axs[1].hist(data['$p_{i}(t_{2})$'], color='k', bins=n_bins)\n",
    "axs[2].hist(data['$p_{i}(t_{3})$'], color='k', bins=n_bins)\n",
    "axs[0].set_xlabel('$p_{i}(t_{1})$', fontsize = labelsize)\n",
    "axs[1].set_xlabel('$p_{i}(t_{2})$', fontsize = labelsize)\n",
    "axs[2].set_xlabel('$p_{i}(t_{3})$', fontsize = labelsize)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, '$p_{i}(t_{3}) - p_{i}(t_{2})$')"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAFgCAYAAACmDI9oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3df7TldV3v8ecLBn9AgJOMUYaOwlUWGrRW0w/LJIRU0EIRwn7c0qBZ1jXu6qqV5s0RswKv4FpZCUmYZrIUycJEgvEypqY2XAUDxsRES/wx5CDhIAi87x/7e/LjZs85333O2fvsc+b5WGuvffbn+/189+fzYfM5r/mez/5+U1VIkiRJGthnpRsgSZIkzRIDsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1Jj5gJzkiCQXJLk+yX1JrlnCsR7eHetLSe5KsiPJLy5jcyVJkrTKrVvpBvTwBOAk4CPAfos9SJKDgA8AdwK/DtwGHAU8aBnaKEmSpDUis36jkCT7VNX93c+XAodU1U8s4jh/CJwKfF9V3bW8rZQkSdJaMfNLLObC8XySPCTJuUn+LcndSa5LctLQbi8ALjIcS5IkaT4zH5B7uhR4PvD7wE8B/wT8bZLvB0jyGOARwO1J3pvkniQ7k5yXxCUWkiRJ+i+rYQ3yvJIcDzwT+Imq2tYV/32SxwG/A5wGHNqVnwtcAjwDOIZBoL4X+M2pNlqSJEkza9UHZOAE4EvAh5K0/dnK4KwyQLrnG6rqV7qf35/kQODlSbZU1e6ptFaSJEkzbS0E5EMYnCH+5oht93XPu7rn/zu0/f3Aq4DDgU9OpHWSJElaVdZCQP4q8AXg2fPs8xngHr51JnnO3OsFvwgoSZKkvcNaCMhbgRcDd1bVjlE7VNU9Sa4CjhvadDywG7h5sk2UJEnSajHzATnJ/gxuFALwSOCgJKd2r98LXAVcCVyV5BzgBuAg4PuBh1TVy7p9zwY+mORi4O3A0cBvA6+uqrun0hlJkiTNvNVwo5CNwGf3sPkxVXVLkgcDLwd+HngUg2UXnwD+qKr+rjnW04E/YHB3vq8AFwKv6XOtZUmSJO0dZj4gS5IkSdO0Vm4UIkmSJC2LmV2DfMghh9TGjRtXuhmSNBHXXnvtbVW1YSXe2/lV0lq2HPPrzAbkjRs3sn379pVuhiRNRJLPrdR7O79KWsuWY351iYUkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktToFZCTHJHkgiTXJ7kvyTXjvEmSfZJsT1JJnrWolkqSJElT0Pcyb08ATgI+Auy3iPc5E/jeRdSTJEmSpqrvEovLq+qwqjoNuGGcN0iyHngN8DvjNk6SJEmatl4BuaruX8J7vBr4ELB1CceQJEmSpmKid9JLcjTwy8DRk3wfSZIkablM+ioWfwS8oapunvD7SJIkSctiYmeQkzwPeDzwU2PU2QxsBnjUox41oZZJ0t7H+VWS+pvIGeQk+wGvBc4B9knyMOCgbvMBSQ4cVa+qLqyqTVW1acOGDZNomiTtlZxftZKSLOohrZRJLbE4gMFl3c4DdnWP67ptlwAfn9D7SpIkSUsyqSUWdwLHDZUdCrwdeDnw/gm9ryRJkrQkvQJykv0Z3CgE4JHAQUlO7V6/t6p2J7kZ2FZVZ1TVvcA1Q8fY2P34yar66FIbLkmSJE1C3zPIjwDeOVQ29/oxwC3dsfZdnmZJkiRJK6NXQK6qW4B5V8tX1calHkOSJElaaZO+DrIkSZK0qhiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJavQKyEmOSHJBkuuT3Jfkmh51fjDJxUluTrI7yaeSvDLJQ5bcakmSJGlC1vXc7wnAScBHgP161jkdOBw4B/g0cDTw6u75ueM1U5IkSZqOvgH58qr6G4AklwKH9Kjzh1V1W/P6miTfAC5I8uiq+tyYbZUkSZImrtcSi6q6f9wDD4XjOR/vnr9n3ONJkiRJ0zDtL+k9Cbgf+MyU31eSJEnqZWoBOcmhwCuAt1bVV6b1vpIkSdI4phKQkzwIeAdwJ/Ab8+y3Ocn2JNt37tw5jaZJ0l7B+VWS+pt4QE4S4C10V8Koql172reqLqyqTVW1acOGDZNumiTtNZxfJam/vlexWIrXAycDP1lVO6bwfpIkSdKiTTQgJ3kZ8CLgZ6rqg5N8L0mSJGk59ArISfZncKMQgEcCByU5tXv93qraneRmYFtVndHV+Tng94E3A19I8iPNIT9TVS6CkyRJ0szpewb5EcA7h8rmXj8GuKU71r7N9qd1z8/vHq0XMAjOkiRJ0kzpFZCr6hYgC+yzcej183lgMJYkSZJm2rRvFCJJkiTNNAOyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUWLfSDZBmRZJF1auqZW6JJElaSZ5BliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqdErICc5IskFSa5Pcl+Sa3rWOzjJxUl2JflakrclefiSWixJkiRNUN/rID8BOAn4CLDfGMd/B/A44EzgfuAc4N3Aj49xDEmSJGlq+gbky6vqbwCSXAocslCFJE8CngYcW1Uf6Mq+AHw0yQlVdfUi2yxJkiRNTK8lFlV1/yKOfSLw5blw3B3nY8Bnu22SJEnSzJnkl/SOBHaMKL+p2yZJkiTNnEkG5PXA7SPKd3XbJEmSpJkzU5d5S7I5yfYk23fu3LnSzZGkNcP5VZL6m2RA3gUcPKJ8fbftAarqwqraVFWbNmzYMMGmSdLexflVkvqbZEDewei1xntamyxJkiStuEkG5CuAQ5M8ea4gySbgsd02SZIkaeb0ug5ykv0Z3CgE4JHAQUlO7V6/t6p2J7kZ2FZVZwBU1T8m+XvgLUlewrduFPJBr4EsSZKkWdX3RiGPAN45VDb3+jHALd2x9h3a53TgfODPGZytfg9w1mIaKkmSJE1Dr4BcVbcAWWCfjSPKbgde0D0kSZKkmTdTl3mTJEmSVpoBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpEavgJzkqCRbk+xOcmuSs5Ps26PepiR/n+Sr3ePqJD+89GZLkiRJk7FgQE6yHrgaKOBk4GzgxcCrFqh3WFdvHfDfu8c64Kokj15asyVJkqTJWNdjnxcCDwVOqao7GATcg4AtSc7tykZ5JnAg8Jyq+hpAkg8DtwEnAX+65NZLkiRJy6zPEosTgSuHgvAlDELzsfPU2w+4F/h6U3ZnV5Yx2ylJkiRNRZ+AfCSwoy2oqs8Du7tte/Kubp/XJXlEkkcA5wO7gHcurrmSJEnSZPUJyOuB20eU7+q2jVRVtwLHAc8Fvtw9TgGeXlU7x2+qJEmSNHkTu8xbku9mcKb4WgbLNE7sfv67JI/aQ53NSbYn2b5zpxlakpaL86sk9dcnIO8CDh5Rvr7bticvZbAO+dSqel9VvY/B2eT7gJeMqlBVF1bVpqratGHDhh5NkyT14fwqSf31Ccg7GFpr3F3CbX+G1iYPORK4oaq+OVdQVfcANwCHj99USZIkafL6BOQrgKcnObApOx24C9g2T73PAU9M8qC5giQPBp4I3DJ+UyVJkqTJ6xOQ3wjcDVyW5IQkm4EtwHntpd+S3Jzkoqbem4DvAf46yTOTPAt4N/DdwIXL1QFJkiRpOS0YkKtqF3A8sC9wOYM76J0PvHJo13XdPnP1rgWeweBmIW8F3sJgWcZPVtV1y9F4SZIkabn1uZMeVXUj8NQF9tk4omwrsHVRLZMkSZJWwMQu8yZJkiStRgZkSZIkqWFAliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqWFAliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqWFAliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqWFAliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqdErICc5KsnWJLuT3Jrk7CT79qx7SpJ/SnJXkv9I8r4kByyt2ZIkSdJkLBiQk6wHrgYKOBk4G3gx8Koedc8E/gq4AjgROBP4NLBu8U2WJEmSJqdPUH0h8FDglKq6A7gqyUHAliTndmUPkOQQ4Hzg16vqz5pNf73URkuSJEmT0meJxYnAlUNB+BIGofnYeer9TPf8F4tsmyRJkjR1fQLykcCOtqCqPg/s7rbtyQ8DnwLOSPLvSb6Z5KNJfnTRrZUkSZImrE9AXg/cPqJ8V7dtTw4FHg+8Avgt4KeArwPvS/Jdoyok2Zxke5LtO3fu7NE0SVIfzq+S1N8kL/MW4DuAM6rqbVX1PuDZwH3Ai0ZVqKoLq2pTVW3asGHDBJsmSXsX51dJ6q9PQN4FHDyifH23bb56BVwzV9CtY74WOKp/EyVJkqTp6ROQdzC01jjJYcD+DK1NHnITg7PIGSoPcP8YbZQkSZKmpk9AvgJ4epIDm7LTgbuAbfPUe0/3fNxcQZKDgR8ArhuznZIkSdJU9AnIbwTuBi5LckKSzcAW4Lz20m9Jbk5y0dzrqtoO/A1wUZJfSvJM4G+BbwJ/vIx9kCRJkpbNggG5qnYBxwP7ApczuIPe+cArh3Zd1+3T+gXg3cB5wKUMwvFTu2NKkiRJM6fXLZ+r6kbgqQvss3FE2Z3Ar3YPSZIkaeZN8jJvkiRJ0qpjQJYkSZIaBmRJkiSpYUCWJEmSGr2+pCetJsnwvWkkSZL68wyyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSY1eATnJUUm2Jtmd5NYkZyfZt++bJNknyfYkleRZi2+uJEmSNFnrFtohyXrgauBG4GTgcOB1DML1K3q+z5nA9y6yjZIkSdLU9DmD/ELgocApVXVVVb0ReBXwv5IctFDlLmC/BvidJbVUkiRJmoI+AflE4MqquqMpu4RBaD62R/1XAx8Cto7fPEmSJGm6+gTkI4EdbUFVfR7Y3W3boyRHA78MvGSxDZQkSZKmqU9AXg/cPqJ8V7dtPn8EvKGqbu7TmCSbuy/zbd+5c2efKpKkHpxfJam/iV3mLcnzgMcDv9e3TlVdWFWbqmrThg0bJtU0SdrrOL9KUn99AvIu4OAR5eu7bQ+QZD/gtcA5wD5JHgbMfaHvgCQHLqKtkiRJ0sT1Ccg7GFprnOQwYH+G1iY3DmBwWbfzGIToXcB13bZLgI8vprGSJEnSpC14HWTgCuClSQ6sqv/syk4H7gK27aHOncBxQ2WHAm8HXg68fxFtlSRJkiauT0B+I3AWcFmSc4DHAluA89pLvyW5GdhWVWdU1b3ANe1BkmzsfvxkVX10yS2XJEmSJmDBgFxVu5IcD7wBuJzBFS3OZxCSh4/V+/bTkiRJ0izqcwaZqroReOoC+2xcYPstQPo2TJIkSVoJE7vMmyRJkrQaGZAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpIYBWZIkSWoYkCVJkqSGAVmSJElqGJAlSZKkhgFZkiRJahiQJUmSpMa6lW6AJEmafUlWugnS1HgGWZIkSWoYkCVJkqRGr4Cc5KgkW5PsTnJrkrOT7LtAnR9McnGSm7t6n0ryyiQPWZ6mS5IkSctvwTXISdYDVwM3AicDhwOvYxCuXzFP1dO7fc8BPg0cDby6e37uklotSZIkTUifL+m9EHgocEpV3QFcleQgYEuSc7uyUf6wqm5rXl+T5BvABUkeXVWfW1rTJUmSpOXXZ4nFicCVQ0H4Egah+dg9VRoKx3M+3j1/T+8WSpIkSVPUJyAfCexoC6rq88Dubts4ngTcD3xmzHqSJEnSVPQJyOuB20eU7+q29ZLkUAZrlt9aVV/Zwz6bk2xPsn3nzp19Dy1JWoDzqyT1N5XLvCV5EPAO4E7gN/a0X1VdWFWbqmrThg0bptE0SdorOL9KUn99vqS3Czh4RPn6btu8Mrj1zluAJwA/VlUL1pEkSZJWSp+AvIOhtcZJDgP2Z2ht8h68nsHl4X6yqvrsL0mSJK2YPkssrgCenuTApux04C5g23wVk7wMeBHwC1X1wUW3UpIkSZqSPgH5jcDdwGVJTkiyGdgCnNde+q27Y95FzeufA36fwfKKLyT5kebhAjhJkiTNpAWXWFTVriTHA28ALmdwRYvzGYTk4WO1t59+Wvf8/O7RegHw5nEbK0mSJE1anzXIVNWNwFMX2Gfj0Ovn88BgLEmSJM20qVzmTZIkSVotDMiSJElSw4AsSZIkNQzIkiRJUsOALEmSJDUMyJIkSVLDgCxJkiQ1DMiSJElSw4AsSZIkNQzIkiRJUsOALEmSJDUMyJIkSVLDgCxJkiQ1DMiSJElSw4AsSZIkNQzIkiRJUsOALEmSJDUMyJIkSVLDgCxJkiQ1DMiSJElSw4AsSZIkNQzIkiRJUsOALEmSJDUMyJIkSVLDgCxJkiQ1egXkJEcl2Zpkd5Jbk5ydZN8e9Q5OcnGSXUm+luRtSR6+9GZrrUuy6IckSdJSrFtohyTrgauBG4GTgcOB1zEI169YoPo7gMcBZwL3A+cA7wZ+fPFNlmbLYkN5VS1zSyRJ0nJYMCADLwQeCpxSVXcAVyU5CNiS5Nyu7AGSPAl4GnBsVX2gK/sC8NEkJ1TV1cvTBUmSJGn59FlicSJw5VAQvoRBaD52gXpfngvHAFX1MeCz3TZJkjRlq2n52mpqq9aWPmeQjwTe3xZU1eeT7O62XT5PvR0jym/qtmkVccJZfi7NkLQUzst75vyqpeoTkNcDt48o39VtW0y9x46qkGQzsLl7eWeST43Y7RDgtnneV9/iWI1nVYzXDP1SXBXjNSNGjdWjp9mAHvOr/z3H55iNZ+bHa4bm1zkzP2Yz6BCWYX7tE5CnpqouBC6cb58k26tq05SatKo5VuNxvMbjePU3C2O10Pw6C21cbRyz8The43PMxteN2calHqfPGuRdwMEjytd325a7niRJkrRi+gTkHQytGU5yGLA/o9cY77FeZ09rkyVJkqQV1ycgXwE8PcmBTdnpwF3AtgXqHZrkyXMFSTYxWH98xSLaOmfeJRj6No7VeByv8The/a2GsVoNbZw1jtl4HK/xOWbjW5Yxy0Lf2OxuFHIj8M8MbvTxWOA84PVV9Ypmv5uBbVV1RlN2JfDfgJfwrRuFfKWqvFGIJEmSZtKCZ5CrahdwPLAvg0u6vQo4H3jl0K7run1apzM4y/znwFuAa4HnLK3JkiRJ0uQseAZZkiRJ2pv0WYO8YpKcnuSyJF9MUkmeP0bdH0vy0STfSPLZJGdNsKkzI8mvJPl01+9rkxzfo86WbnyHH8+YRpsnLclRSbYm2Z3k1iRnJxn+a8eoegcnuTjJriRfS/K2JA+fRptX0mLGK8nGPXyGLplWu1dCkiOSXJDk+iT3JbmmZ70V/Ww5ty6O8+tozrHjc54d37Tn25m6DvIIpwIbgfcAZ/atlOQI4Mqu3suAHwLOS7K7qt40gXbOhCQ/C7wR2AJ8EHgB8J4kP1hV/7xA9a8BwxP2TcveyCnr1tBfzWAd/cnA4cDrGPzj8BXzVAV4B/A4Bp+9uTX07wbW7Br6JY4XDL5v8KHm9Vq/wP0TgJOAjwD7jVFvpT9bzq1jcn4dzTl2fM6zizbd+baqZvYB7NM9fwdQwPN71rsA+BdgXVP2J8C/0S0rWYsP4FPAn7fjB3wS+MsF6m0Bblvp9k9oTF7G4LrbBzVlvwnsbstG1HtS95l7SlP2Q13ZCSvdrxkcr43d2Dxrpfsw5fHap/n5UuCaHnVW/LPl3LqoMXN+Hd0/59jpjdleOc82/Z/qfDvTSyyq6v5FVj0RuKyq7m3KLgG+F3jikhs2g5I8lsG/kN4xV9aN3zsZjMfe6kTgyqq6oym7BHgocOwC9b5cVR+YK6iqjwGfZW2P52LHa6+0yDlqxT9bzq3jcX6dl3Ps+JxnF2Ha8+1MB+TFSHIAcBgPvBnJ3J+zRt28ZC2Y69eofn9nkg0L1H9YktuSfDPJx5OcsvxNXBEPuDFNVX2ewb/U5/ss7OmGNjctUG+1W+x4zbm4Wxv2xSTnJXnoJBq5yq3Kz9ZePLeC8+t8nGPH5zw7PYv+nK25gAw8rHu+fah87vbW66fYlmma69di+n0zgz/vnAY8F7gVeNcamcTX88AxgcG4zDcmi6232i2233cDfwycweCykBcAv8rgrIi+3Wr9bO2tcys4v87HOXZ8zrPTs+jP2VS/pJfkYOC7F9qvqrwVNdMbr6r6y6H3vRz4MPC7wGVLObb2DlX1ReBFTdE1Sb4M/EmSY6rquhVq2l7BuXV8zq9abZxnp2vaV7E4DfizHvtlCe8x9y+Fg4fK5/6lsIvVY5zxmuvXwXz7v5bG7ndVVZLLgHOS7FtV9/WtO4N28cDPAgzGZb4x2QWM+rPpQvVWu8WO1yiXMvgC1w8ATtzfMonPlnPr+Jxfl4dz7PicZ6dn0Z+zqS6xqKo3VVUWeizxPb7O4BvVw2tL9rSGbGaNOV5z/RrV769W1c5x3757rHY7GBqTJIcB+zP/Z+EB9Tp7Ws+0Vix2vEapoWcNLPtny7l1fM6vy8Y5dnzOs9Oz6M/ZWlyDDHAF8Jyhi26fzmByX+h6latSVf0rg8svnTZXlmSf7vUV4xwrSRislbtuDZzduAJ4epIDm7LTgbsY3AZ9vnqHJnnyXEGSTcBjGXM8V5nFjtcop3bP1y5Hw9aQ1fzZ2uvmVnB+XYBz7PicZ6dn8Z+zSVyrbrkewFEM/uP/AoN/Hb2he31ss8+xwL1DZUcAdwJ/BRzH4AsS3wTOXOk+TXi8fha4j8GFxo8D3szgf7gnLjBe24CzgKcBzwHey+Bi2j+90n1ahjFZD3wRuAo4AdjcfTZ+b2i/m4GLhsquBP4VOAV4NoProP7DSvdpFseLwbVeX9eN1QnA2d1n710r3acJj9f+3Zx0KvCPwA3N6/1n9bPl3LqoMXN+HT0uzrFTGrO9dZ5t+j/V+XbFO7zAYGzhW3+Kah/XNPv8RFf2E0N1nwx8DPgGcAtw1kr3Z0pj9ivdB+Ru4P8Bxw9tf8B4ARd1H567gK8D/wCcuNJ9WcYxOQp4f9e/LwKvBvYd2ucW4M1DZQ8DLmaw5vAOBqHgkJXuzyyOF/A8YDuDO4bd030GzwYevNL9mfBYbdzDHFXAxln9bDm3LnrcnF9Hj4tz7BTGbG+dZ5v+T3W+TVdZkiRJEmt3DbIkSZK0KAZkSZIkqWFAliRJkhoGZEmSJKlhQJYkSZIaBmRJkiSpYUCWJEmSGgZkSZIkqWFA1szIwCeS/NJQ+eOSbEnyoJ7HeUOSiybTysVb6/2TNJvW+tyz1vunlWFA1iz5GeA7GdwGsnUi8GtVdU/P4/wf4OeTHLGcjVsGa71/kmbTWp971nr/tAIMyJolZwFvrapvDpUfDfxz34NU1S3AB4FfXb6mLYu13j9Js2mtzz1rvX9aAQZkTUySryZ5dZK3JfmP7vUFSR48Yt8jgB8FLh0q/wrwy8BxSap7PLXH27+LwZmAiX3G13r/JM2mtT73rPX+aXXwA6CJSPIoYD3wG8BXgdOA1wObgf8xosrxwNeB65pjBPhp4BvAa4AndY8P92jCh4HvAr5v0Z2Yx1rvn6TZtNbnnrXeP60e61a6AVqzjume31BVv939/P4kTwGeAZw3tP8PADdV1f1zBVVVSb4EPAR4T1V9ZIz3vwG4D/ghmolzGa31/kmaTWt97lnr/dMq4RlkTcoxwG7gtUPl/8LgyxQkOTLJm7ryQ4HbRhznaKCAT84VJDksydYkNyW5Icm53RmD/1JV9wK3d8d9gCQHd+8/72Ml+tfV/dMkX0hSo958of5JWrMmPfdsS3JdkuuTXJrkoHb7ap5bl+N3h/YenkHWpBwDfKCq/mOo/FDgVoCq2gGc2ZU/hMGkOOxo4F+r6utN2b3Ab1XV9gwu33MVcAqDtWOtu7vjjnIa8Gc9+pE9lE+yfwBvB7YAX5qnbfP1T9LaNOm556er6msASc4DXgr876F9Vuvcuhy/O7SX8AyyJuUYuslsTpIDgacCW7vXb0zyvG7zV4GHjTjO0Qz9mauqvlhV27uf7wGuBw4bUfdh3XEfoKreVFVZ6LES/eva94Gq+vI87z9v/yStWZOee+bC8T7AAQzOwg5blXPrcvzu0N7DgKxll+QA4PDu0Tqre/6L7vn7gU90P38KeMyIwx0J7JjnvR4OPBu4cqh8A7A/gz/LLatp9m+eNkysf5Jm07TmniTvBb4MPB44d2jbmphbV+J3h1YXA7Im4fsYfLYOT/K7SY5L8gfAK4HNVXV7d3bisXxrEvoQ8KhucmrdATwlyVOS/Ei7Xqy75M+lwOur6qahepsYnPno863lcU2lfwuYZP8kzaapzD1VdRKDJQ0fA35tqN6qn1tX8HeHVhEDsibhGAaT0zOBk4H3dT+fVgoI15IAAAF4SURBVFXv6PZ5HIP1YXPfPL6GwZ+0njF0rN9m8MWMrcC7q6oAkuwLvA34eFW9bkQbngFsG7GObTlMvH89TLJ/kmbT1OaeqrqPwRnbXxyqt6rn1hX+3aFVxC/paRKOAT5ZVdczuATPKO2fyKiqe5L8JfA84K1N+QeBJ4yofwHwn8CLhzd0E+BzGUyQkzCN/u3RFPonaTZNdO5Jsh54UPP9h+fS3IlujcytK/m7Q6uIZ5A1Ccew8PUjv22S67yWwV2PHjdfxSQ/BpzB4E9hH0/yiSRnNbucBtwFXDJWq/ubaP8Akrwpyb93P/97c0kjmHz/JM2mSc8964H3dJd4+yTwROB/NttX9dw6A787tIqk/190pYV167y+Bry0qi5YRP3nAV+sqm1LaMPPAl+oqg8s9hjzHHtN90/SbFrrc89a759WHwOyJEmS1HCJhSRJktQwIEuSJEkNA7IkSZLUMCBLkiRJDQOyJEmS1DAgS5IkSQ0DsiRJktT4/zTiXT5YWsc5AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x360 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axs = plt.subplots(1, 2, sharey=True, tight_layout=True, figsize=(10, 5))\n",
    "labelsize = 15\n",
    "plt.rcParams['xtick.labelsize'] = labelsize\n",
    "plt.rcParams['ytick.labelsize'] = labelsize\n",
    "plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))\n",
    "n_bins = 20\n",
    "axs[0].hist(data['$p_{i}(t_{2})$'] - data['$p_{i}(t_{1})$'], color='k', bins=n_bins)\n",
    "axs[1].hist(data['$p_{i}(t_{3})$'] - data['$p_{i}(t_{2})$'], color='k', bins=n_bins)\n",
    "axs[0].set_xlabel('$p_{i}(t_{2}) - p_{i}(t_{1})$', fontsize = labelsize)\n",
    "axs[1].set_xlabel('$p_{i}(t_{3}) - p_{i}(t_{2})$', fontsize = labelsize)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Vertex degree distribution "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "count   1660927.000\n",
       "mean         17.482\n",
       "std          23.241\n",
       "min           1.000\n",
       "25%           4.000\n",
       "50%          11.000\n",
       "75%          22.000\n",
       "max        1315.000\n",
       "Name: fr_number, dtype: float64"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data['fr_number'].describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Histogram (both axes in logarithmic scale)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWIAAAFJCAYAAAChNyePAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3df5xddX3n8debSJBQgqMkdMFkg8QSEhZUxkUWEQUeKpEIxVBWdLUWiVRb2IrIo91UElrchS2hBcKPCCIiFSSG6lhYCmgiFVESwAAhLa6GbApIpAOJBoSEz/5xvoM3lzuZe2bOme+dmffz8bgP55zzvfd+5jziZ758zveHIgIzM8tnp9wBmJmNdU7EZmaZORGbmWXmRGxmlpkTsZlZZq/JHUCn2XPPPWPatGm5wzCzUWbVqlW/jIhJra45ESeS5gBzpk+fzsqVK3OHY2ajjKTH+7vm0kQSET0RMW+PPfbIHYqZjTFOxImkOZKWPPfcc7lDMbMxxok4cY/YzHJxIjYzy8yJOHFpwsxycSJOXJows1yciM3MMnMiNjPLzIk4cY3YzHJxIk5cIzazXJyIzcwy81oTFVi4cOGQ3n/uuedWFImZjUTuEZuZZeZEbGaWmRNx4lETZpaLE3HiURNmlosTsZlZZqN61ISkdcAW4MV06pSIWJMvIjOzVxvViTiZHRHrcgexI4MZ/uYhb2ajR8eVJiRNl3SVpNWStkla3k+7mZLukrRF0hOSzpM0bpjDNTMbsk7sEc8CZgP3Aju3aiCpC7gTWAMcD+wHXETxh2V+U/N/kCTgO8CCiHipprjNzAal43rEQE9ETImIk4BH+mlzOrArcGJE3BERVwILgc9KmtjQ7p0R8RbgcGAm8Lk6AzczG4yOS8QR8XIbzY4Fbo+ITQ3nbqRIzkc2fNaG9L+/Aq4B/kuFoZqZVaLjEnGbZgBrG09ExHqKERIzACTt1tc7lvQa4EPA6lYfJmmepJWSVm7cuLHWwM3MmnVijbgdXcCzLc73pmsAewHLJO0EjAN+CJzf6sMiYgmwBKC7uzsqjzajoS5IBB6hYVa3kZqIBxQRPwPe0m57SXOAOdOnT68vKDOzFkZqIu4FWs1F7krXRr0qerpm1hnaTsSSLml1PiLOqC6ctq0l1YL7SJoCTKCpdtyuiOgBerq7u08benhmZu0b8GGdpDPTj8cB7wIeBlY1vHK4DXifpN0bzp0MPA+sGMwHevU1M8ulnVETp6T//T3gS8BpFBMtro+I66oOSNIESXMlzQX2ASb1HUuakJpdCfyG4mHcMZLmAQuARU1D2trm1dfMLJd2EvFmgIjYGhGLKcbpTgLuScmyapOBm9PrHRQTMfqOJ6dYeoGjKUZD9FBM5rgYGPTjffeIzSyXdmrEpwBIegjoG9oliodlN1Ekw8qkBXrURrs1wFEVfq9rxEPgYXJmgzdgIo6Ip9OPx9Uci3Uoj9Awq1eZmXXrgSOAj0fE4xS9471qiSoDlybMLJcyifhy4DDgw+l4M7C48ogy8cM6M8ulzISOQyPibZIegOKBmaTxNcVlZjZmlEnEL6WF1wNA0iSgnZXSRgRPcc7PO5XYWFWmNHEJcAswWdL5wD8DX6wlqgxcmjCzXNruEUfEDZJWUYzfBTg+IgY1ndjMzH6r7R6xpMOBvwTuBw4ArpB0WF2BmZmNFWVKE5cByylmst1NUaoYNaMmPHzNzHIpk4hfjogvAc9ExE0RcQttzIAbKVwjNrNcyoya2Jb+9w8A0s4XI3WrJRslPNLCRoMyiXQ2QET8JB1PAOZVHpGZ2RhTJhF/WFLffnBExK8i4kc1xGRmNqaUScR7AfdJ+oak90saNfVh8MM6M8un7UQcEfOBNwPXAH8IPCbpi5L2qym2YeWHdWaWS6mHbRERwFPptZVis86lki6sITYzszGhzOahZwIfA34JXA2cHREvpdETjwGfrydEM7PRrczwtdcDJ6a1iF8RES9L6uhF4yUtBj4dEaOqrm1mo0OZtSb6HXwZEY9WE071JB0B/E7uOMzM+lNmrYnXSvqspGWSvinpzyS9tuqAJE2XdJWk1ZK2SVreT7uZku6StEXSE5LOS8t0NrbZBfhfwOeqjtPMrCplShNfpdiV49J0fApwPXBSxTHNopg8ci+wc6sGaTzzncAa4HhgP+Aiij8s8xuafgG4JiI2jrLRdmY2ipRJxAdGxMyG4+9JWlN1QEBPRHwLQNJSYM8WbU4HdqWoWW8C7pA0EVgg6cKI2CTpIOBQtk/MZmYdp0wivl/SOyLiXgBJhwIrqw4oItrZ9eNY4PaUhPvcCFwAHEmxQtzhwEzg5329YUnrgLdHxMYqY7aRZbC7UnuNCqtLmXHEhwD3SFqXEtoPgbdLekjS6lqi698MYLtF6SNiPbAlXSMiroiIvSNiWkRMS+emtUrCkuZJWilp5caNztFmNrzK9IjfX1sU5XUBz7Y435uulRIRSyQ9CcwZP378IUMNzsysjDJTnB9vfgH/iaIU0PEGGkPsKc5mlkuZmXXfbj4FvBP4CPCbKoNqQy/QKmN2pWuleRdnM8ulTGniAOCTDccCZkTErdWG1Ja1pFrwK8FIUyjWSB7UhqYR0QP0dHd3nzb08MzM2lcmEW+OiBWNJyRtrjiedt0GnC1p94joi+Fk4HlgRf9v6597xGaWS5lRE7Mk/VTSj9Psuj8C6phZN0HSXElzgX2ASX3HkiakZldSlEOWSTpG0jxgAbCoaUhb21wjNrNcyvSI9wbGUazbsC/FjLr9Jb0LWBMRv6wopsnAzU3n+o73BdZFRK+koyl2lu6hGEFxMUUyHhT3iM0slzKL/jyTfnwa+BlwVxo//B6KpTErScQRsY42doeOiDXAUVV8Z/o814jNLIsyPeJXiYjLqwokN/eIzSyXUjt0jGauEZtZLkPqEZuNJYNZo8LrU1g7yqxH/Nctzo1r1XYk8i7OZpZLmdLEPpI+3HcgaTLFmsCjgksTZpZLmdLEp4DbJf1fIIBrgXNqicrMbAwZMBFL+ipwP/AA8Bng74GtwAkR8dN6wzMzG/3aKU18hWJc7yeArwHTKBbW+Wia/TYquEZsZrkM2COOiO8C3+07lvQaigWADqbYimhpbdENI0/oMLNcSg9fi4itwEPp9bXKIzIzG2M8ocPMLDMnYjOzzMpM6JjZ4ty7K43GzGwMKtMj/oakc1TYVdKlwP+sK7Dh5lETZpZLmUR8KDAFuAe4D3gCOLyOoHLwzDozy6VMIn6JYiuiXSl25vh5RLxcS1RmZmNImUR8H0UifjtwBPBhSc07aXQUSSsk/UTSaklLJU3MHZOZWbMyifjUiPhCRLwUEU9GxPHAt+sKrCIfjIiDI+IgYD1wdu6AzMyalZnQ8TFJH2tx/vqqggGQNJ0iYR4GzALujoh3t2g3E7g0tXsWuBpYGBHb+tpExHOp7U7AbsCvqozVzKwKZRLxB4DNwOUUOyjXZRYwG7gX2LlVA0ldFEtwrgGOB/YDLqLo4c9vansrRTnlEeCs2qI2a8GLyVs7ypQm9ge+BJxGkSCvj4jraoipJyKmRMRJFMmzldMpHhqeGBF3RMSVwELgs8114IiYDfwu8GPg0zXEa2Y2JG0n4ojYGhGLgSOBScA9day+1uZIjGOB2yNiU8O5GymS85EtPnMbcB3QqrRiZpZV26UJSQ9RLAgPxbKYewA3ATm2S5pBw4pwABGxXtKWdK0nlS/GR8QvUpMPAQ+3+jBJ84B5AFOnTq0taDOzVsrUiI+rLYryuige0DXrTdf62twkaReKPxyPAn/a6sMiYgmwBKC7uztatTEzq0vbiTgiHq8zkKpFxM8oHtK1RdIcYM706dPrC8rMrIUypYl3tTh9KfAMcE1E3FBZVAPrpSiNNOtK18zMRowypYlWkyH2iYiDqwqmhLUUteBXSJoCTEjXSvMOHWaWS5nSxJzmc5Lurjactt0GnC1p94jYnM6dTDEFe8VgPtClCTPLZagLw1f+YEvSBElz09C4fYBJfceSJqRmV1JMKlkm6Zg06mEBsKhpSFvbvPqameVSpka8me0TryhWYavaZKB5MaG+432BdRHRK+lo4DKgh2IExcUUyXhQ3CM2s1zKlCZ2bz5XR2kiItZRJPmB2q0Bjqrwe10jNrMsOq40YWY21gxmZp0a/ndaPWENP5cmzCyXkTqzrnIuTZhZLmUS8fMUS00ekY6X44XWzcyGrEyN+FqK7ZIOSq9VwJfrCCoH7+JsZrmUScRvjIhLImJTel1KsavzqOBxxGaWS5lEvFXSB/oOJM2m2LHDzMyGoEyN+NPA30v6KrCNYm3fP6olKrMxzNsrjT1lJnT8CNhP0oPA29rcSWPE8PA1M8ul7dKEpM2SNgEHAs82HI8KrhGbWS5DmuJsZmZDV6ZHLEkflfSX6XiKpP9cX2hmZmNDmVETlwOHAaek418BiyuPyMxsjCmTiA+NiM8ALwBERC8wvpaoMvCEDjPLpUwifknSONKKa5ImAaNm5IQf1plZLmUS8SXALcBkSecD/wx8sZaozMzGkDKjJm6QtAo4mmIJzBMi4tHaIjMzGyPKzKwjItYyyF2Sh1va1fkrwN4UJZR/BM6JCC9mb2YdpczC8C1XWouITp3mvJUi8a6UNB64AzgR+GbesMzMtlemRvw+YAbFUpj/2PCqlKTpkq6StFrSNknL+2k3U9JdkrZIekLSeelhIgAR8WRErEw/vwisZhStFmdmo0eZ0sQU4P3AfwPGAddGxG01xDQLmA3cC+zcqoGkLuBOYA1wPLAfxaL1OwHzW7R/A3AC8N4a4jUzG5K2e8QR8XJE3Ar8FbAF+JOaYuqJiCkRcRLwSD9tTgd2BU6MiDsi4kpgIfBZSRMbG0raBVgK/K0fLppZJyozxXmepFspkuDfRcQHBnrPYLS5qtuxwO0R0bjo0I0UyfnIvhOpVHED8EBEXFRpoGZmFSlTI74SmA68B7gu1XBX1xPWgGbQNHojItZT9NRnNJy+imLx+rN29GHpj8xKSSs3btxYdaxmZjtUpka8b21RlNcFPNvifG+6hqTDgVMpFrB/QBLAlyPikuY3RcQSSU8Cc8aPH39IbVGbmbVQJhH316s8o4pAqhYRP6CYeNJu+x6gp7u7+7T6ojIze7Uyifg4YBPFKmy/qSectvUCrRaF6ErXSvMOHWaWS5ka8e8BXwJOoxhWdn1EXFdLVANby/a14L6ZdBMY5Mw/L/pjZrmUWWtiK7BY0rXAmcA9kv4mIpbWFl3/bgPOlrR7RPTtJH0y8DywYjAf6B6xjWTecHRkKzPF+SHSEpgUtdc9gJsoJndURtIEigkdAPsAEyXNTce3RsQWihEcZwDLJF0AvAlYACxqGtLWNteIzSyXsjXi4TAZuLnpXN/xvsC6iOiVdDRwGdBDMYLiYopkPCjuEZtZLmVKE4/XGUjD96yjjdEOEbEGOKrC73WP2MyyKPOwblTzVklmlosTceJRE2aWixOxmVlmg07Ekv5DWtlsVHBpwsxyGUqP+HpgraS/qSqYnFyaMLNcSu1Z1ygijlGxks7MCuMxMxtzykzoeNWqZQAR0ZGL/piZjRQDliYknZl+PA54F8WykqsaXqOCa8Rmlks7NeJT0v920qI/lXON2MxyaScRb4Zi0Z+IWEyxFdEkikV/5u7wnWZmNqB2asSnwPAt+mNmNtYMmIgj4un043At+mNmNqaUGUe8HjgC+HhaACiAvWqJysxsDCmTiC8HDgM+nI43A4srjygTj5ows1zKJOJDI+IzwAsAEdELjK8lqgw8asLMcimTiF+SNI70wE7SJODlWqIyMxtDykxxvgS4BZgs6XxgLjC/lqgqIukK4IPA3hEx4GLzZmOJ97nrHGV26LhB0irgaIrhaydExKO1RVaNr1Nsn/RU5jjMzPpVavW1iFgbEYsj4rK6krCk6ZKukrRa0jZJy/tpN1PSXZK2SHpC0nmpdNIY7/cj4hd1xGlmVpVOXPRnFsUuzvdSTKVuFUsXcCewBjge2A+4iOIPS0eXS8zMmpWpER8PfKGuQBr0RMS3ACQtBfZs0eZ0YFfgxIjYBNwhaSKwQNKF6ZyZ2YhQJhE/MxyL/EREOyMxjgVub0q4NwIXUKyF0VNHbGZmdShTI54h6UFJ90paJuksSa+tLbIBYgHWNp6IiPXAlnStFEnzJK2UtHLjxo0VhWhm1p4yifgA4ETgVOAKYH/g6jqCakMX8GyL873pGgCSrpa0If28QVLLeCNiSUR0R0T3pEmTagnYzKw/ZYavPd5w+AhFXfaC6kOqTkR8st22kuYAc6ZPn15jRGZmr9Z2j1iFj0r6QjqeCnyztsh2rJdiGc5mXemamdmIMVIX/VlLUy1Y0hRgAk2143Z5rQkzy2WkLvpzG/A+Sbs3nDsZeB5YMZgP9OprZpZLxy36I2mCpLlpG6Z9gEl9x5ImpGZXAr8Blkk6RtI8iqnMiwY7htg9YjPLpRMX/ZkM3Nx0ru94X2BdRPRKOhq4jGLM8LPAxRTJeFD8sM7Mcum4RX8iYl36/IHarQGOqvB7e4Ce7u7u06r6TDOzdpRZa2JBRCxgkA/DzMystTI14g/WFkUH8MM6M8ulTI14sqTPNp+MiEUVxpONSxNmlkuZRDwO+B3aqN+amVn7yiTipyLivNoiycyjJswslzI14jtqi6IDeByxmeVSZvja59POGG8GXttw/vt1BGZmNlaUGb72SeBM4I3Ag8A7gB9S4VheM7OxqEyN+Ezg7cC9EfEeSTOAL9YT1vBzjdhsYAsXLhzU+84999yKIxldytSIX4iIFwAk7RIRaykWhx8VXCM2s1zK9Ig3SHod8A8Ui8L3Ao8P8B4zMxtAmYd1v59+XCDpexQLs/+fWqIyMxtDyvSIG/0U6IqIF6sMxsxsLCqzVdL/lvS0pP8B/BNwg6SL6wtteHmtCTPLpczDut8HDgQ+BxwCvA14bx1B5eCHdWaWS5lEvCkinqZYmP2FiNhGsUuGmZkNQZlEPEPSamB/SaslPUSHD1+TdKCk+yU9JunbTXvcmZl1hDIP6w6oLYr6XAnMj4hbJV0IfB74y8wxmZltp+0ecUQ8DrwOmJNer0vnKiVpuqSrUq97m6Tl/bSbKekuSVskPSHpvLS5ad/1vYB9I+LWdOoa4ENVx2tmNlRlRk2cCdxAsbnnZOBrkv60hphmAbOBfwH+tZ9YuoA7KXaUPh44DzgLaJx/+UZgQ8PxemBKDfGamQ1JmdLEqcChEfFrAEkXUCz6c2nFMfVExLfSdywF9mzR5nRgV+DEiNhEMdNvIsVkkwvTOS9gb2YjQpmHdQK2NRxvo4ZkFxEvt9HsWOD2lHD73EiRnI9MxxsoesV9prJ9D9nMrCOUScTXAj+StEDSQuBHFHXXHGbQtJt0RKwHtqRrRMRTwDpJs1OTU4FlwxmkmVk7yqw1sSg9OHsnRW324xHxYF2BDaALeLbF+d50rc8fA9dJ+juKmvNHWn2YpHnAPICpU6dWG6mZ2QDKPKw7nOKB2H0U44cvlnRYXYFVISJWR8RbI+LNEXFcRLScvxwRSyge9N0/fvz44Q3SzMa8MqWJy4DlQA9wN3AJsLiGmNrRS7H6W7OudK00T3E2s1zKJOKXI+JLwDMRcVNE3EK+kQlrSbXgPpKmABNoqh23y4v+mFkuZRJx34iJPwCQtFPJ91fpNuB9TVOWTwaeB1YM5gPdIzazXMqMI54NEBE/SccTSA+4qiRpQt93AfsAEyXNTce3RsQWiqnLZwDL0njmNwELgEVNQ9rKfK/3rDOzLBQR5d8kvRnYOSLWVB6QNA34eT+X942IdandTIq69WEUIyiuBhakVeEGrbu7O1auXFnqPYPdUNHM+jfaNhyVtCoiultdK71Dh6S/AD4O/FrSioj4s6EG2Cgl2gFrz+mPwFFVfa97xGaWy2BqvHOBt1AsDv/OasPJxzViM8tlUA/bIuL5KGoaz1ccTzYeNWFmuZSZ0PFQWhh+RsPC8G+tL7Th5R6xmeVSpkZ8XG1RmJmNYWUS8Vn9nD+jikBy88M6M8ulTI34eGBVi9eo4NKEmeVSpkf8TERcV1skZmZjVJlEPEPSg8ALwBPAD4DFEfFCLZGZmY0RZXdxHkexC8bewEkUs9k+WkNcw841YjPLpdQuzhHxs4h4JCLuiIh5wL/VGNuwco3YzHIpNcVZ0sHAEenw7og4p/qQzMzGljITOs4EbgAmp9fXJP1pXYGZmY0VZXrEpwKHRsSvAdLykz8ELq0jMDOzsaLMOGLx28XhST/n2qHDzGzUKNMjvhb4kaRb0vEJwDXVh5SHR02YWS5lRk0sAj4B/Ht6fSIi/rauwIabR02YWS6lRk1ExP3A/TXFYmY2JrWdiCV9udX5iPij6sKplqQrgA8Ce0eE69lm1pHK9IjfDZxdUxx1+TrFpqJPZY7DzKxfZRLxcxHxzaF8maTpFMn8MGAWxaSQd7doN5NiWFzjxqALy24MGhHfT583lLDNzGpVJhGX3+751WYBs4F7gZ1bNZDUBdwJrKFYenM/4CKKB4vzK4jBzKyjlF19bXXDsYCIiINKfEZPRHwLQNJSYM8WbU6nWFjoxIjYBNwhaSKwQNKF6RyS7gemtnj/7RHxkRIxmZllVXb1tSGJiJfbaHYsRTLd1HDuRuAC4EigJ33W24Yaj5lZJ2g7EUfE43UG0mAG8N2m714vaUu61lP1F0qaB8wDmDq1VSfbzKw+ZaY4D5cuigd0zXrTtbZJulrShvTzBklXt2oXEUsiojsiuidNmlQ6YDOzoSg1oWOkiYhPttvWU5zNLJehrkf8k+pDohdoNc+4K10zMxtVOnE94rUUteDG754CTEjXauG1Jswsl05cj/g24GxJu0fE5nTuZOB5YEXF3/UKlybMLJdhXY9Y0gRJcyXNBfYBJvUdS5qQml0J/AZYJumYNKJhAbCoaUhbpdwjNrNchns94snAzU3n+o73BdZFRK+ko4HLKIaqPQtcTJGMa+MesZnlUmYc8SJJy4F3plOfiIgHynxZRKyjjV50RKwBjirz2UMVET1AT3d392nD+b1mZl6P2MwsszLrEV/S6nxEnFFdOPm4NGHWWRYuXFj6Peeee24NkdSvzMO644FVLV6jgh/WmVkuZUoTz0TEdbVFYmY2RpVdBvNB4AXgCeAHwOKIeKGWyIaZSxNmlkuZ0sQBwIkUEzuuAPan2DljVHBpwsxyGewymI9QLNh+QfUhmZmNLWXWmrhO0usajrsArxlpZjZEZUoTB0XEK+sER0Qv8NbqQ8pD0hxJS5577rncoZjZGFMmEe+UesEASHo9o2g9Y9eIzSyXMon0IuCHkvrWhjgJOL/6kMzMxpYyD+u+KmkV8J506sS0JoSZmQ1B2bUmHqEYMWFmZhXxWhOJJ3SYWS5lHtZ9AHgX8DBea8LMrDJlEvH+wJeA04Cdgeu99oSZ2dC1nYgjYmtELAaOpJjIcU/a8qgjSZoi6S5Jj0p6RNKFkkpt7WRmNhzK1IgfAqLvkGLL+5uAcTXEVYWtwDkRsVLSeOAOirUyvpk3LDOz7ZUZNXFcw8+NCbltkqYDZwOHAbOAuyPi3S3azaTYHfowij3rrgYWRsS25rb9iYgngSfTzy9KWg1MKROvmdlwGDARS/o5v028211K599U4vtmAbOBeynqzK2+rwu4E1hDsRj9fhSTSXYC5pf4rsbPfAPFZqfvHcz7zczq1E6PuLvhZwHf5beTOsrqiYhvAUhaCuzZos3pwK4UE0Y2UazyNhFYIOnCdA5J9wNTW7z/9oj4yCsBS7sAS4G/jYhHBxm3mVltBkzEEfFM47Gkrc3n2hURL7fR7FiKZLqp4dyNwAUUDwp70me9baAPkjQOuAF4ICIuKh+xmVn9ygxfQ9KbKFkXHoQZwNrGExGxHtiSrpVxFbAZOKua0MzMqtdOjbhvtMQuwATgUzXH1EXxgK5Zb7rWFkmHU+wm8jDwQBq59uWIeNUMQUnzgHkAU6e2qnaYmdWnnRpx32iJFyLiF3UGU6WI+AFt9t4jYomkJ4E548ePP6TeyMzMtjdgaSIiHk+v4UrCvRRjlJt1pWu18BRnM8ulVI14mKylqRYsaQpFWWRty3dUwDt0mFkunZiIbwPeJ2n3hnMnA88DK+r6UveIzSyXYd3qSNIEigkdAPsAExvWq7g1IrYAVwJnAMvSLtFvAhYAi5qGtFUdm5fBNLMshnvPucnAzU3n+o73BdZFRK+ko4HLKMYMPwtcTJGMaxMRPUBPd3f3aXV+j5lZs2FNxBGxjjZGMqQtmI6qPaAG7hGbWS6dWCPOwjViM8vFiTjxqAkzy8WJOHGP2MxycSI2M8vMiThxacLMcnEiTlyaMLNchnscsZlZbRYuXFj6Peeee24NkZTjHrGZWWZOxIlrxGaWixNx4hqxmeXiRGxmlpkTsZlZZk7EZmaZORGbmWXmRJx41ISZ5eJEnHjUhJnl4kRsZpaZE7GZWWaKiNwxdBRJG4HHW1zaA2guIDef29HxnsAvKwpzoLiqeM+O2vR3baj3CHyf2jkeTfeov+uDOdfp/5ZeFxGTWl6NCL/aeAFLBjq3o2Ng5XDFVcV7dtSmv2tDvUe+T2Pv31K796OdcyP535JLE+3raePcQMd1GMx3tPOeHbXp71qn3qPBfs9Yu0857lF/1wdzbsT+W3JpYphIWhkR3bnj6HS+TwPzPWrPSLpP7hEPnyW5AxghfJ8G5nvUnhFzn9wjNjPLzD1iM7PMnIjNzDJzIjYzy8yJuANIOlDS/ZIek/RtSbvnjqkTSbpC0r9J8oONfkiaIukuSY9KekTShZKUO65OI2mFpJ9IWi1pqaSJOeNxIu4MVwLzI+LNwFrg85nj6VRfB96WO4gOtxU4JyIOAN4KHAqcmDekjvTBiDg4Ig4C1gNn5wzGiXiQJE2XdFX6i7pN0vJ+2s1MPZQtkp6QdJ6kcQ3X9wL2jYhb06lrgA/V/xsMj6ruE0BEfD8ifjEsgQ+zqu5TRDwZESvTzy8Cq4Epw/JL1Kzif0vPpbY7AbsBWf8r6zU5v3yEmwXMBu4Fdm7VQFIXcCewBjge2A+4iOIP4PzU7I3Ahoa3rWeU/B8nqfK21NYAAAW2SURBVOo+jXaV3ydJbwBOAN5bT8jDrtJ7JOlW4O3AI8BZtUXdjjrmYo+FF7BTw89LgeUt2vw50AtMbDj3eWBL3zmgG/hRw/Vdgc25f79Ou09N7SP379Xp9wnYBfgecFbu361T71G6Ng64EPh8zt/NpYlBioiX22h2LHB7RGxqOHcjRbI9Mh1voOgV95nK9j3kEa3C+zSqVXmf0n+G3wA8EBEXVRpoRnX8W4qIbcB1wMcqCXKQnIjrNYPi4dsrImI9xV/nGen4KWCdpNmpyanAsuEMsgMMeJ8MaP8+XQVsJvd/bucx4D2S1JWezfT5EPDwsEXYghNxvbqAZ1uc703X+vwxcL6kx4CZFP+pNJa0dZ8kXS1pQ/p5g6Srhym+TjHgfZJ0OMUf827gAUkPSjpj+ELMrp1/S13Ad9JDv4eAA4Ezhym+lvywrgNExGqKoUa2AxHxydwxdLqI+AHgccM7EBE/o3hI1zHcI65XL8XK/M260jUr+D61x/dpYCPyHjkR12stTTVOSVOACTTVscY436f2+D4NbETeIyfiet0GvK9pyvLJwPPAijwhdSTfp/b4Pg1sRN4j14gHSdIEisHlAPsAEyXNTce3RsQWiqnLZwDLJF0AvAlYACxqGl4zavk+tcf3aWCj+h7lHqQ9Ul/ANIppka1e0xrazQS+S/EX+Ungr4BxueP3feqsl+/T2L5H3qHDzCwz14jNzDJzIjYzy8yJ2MwsMydiM7PMnIjNzDJzIjYzy8yJ2MwsMydiM7PMnIjNzDJzIraOJWmapIcbjveUtC5jSFlI2ivtSnyfpD/rp83Rkq4f7tisGl70x6zDRcQvgKMHaHYw8MAwhGM1cI/YOtkLwPjmk409ZUk7S/qZpMsarn8sbYPzk75eYnrP82nroPVN7T8r6eH0+u+tAkm98RfT+38q6TsNn7tW0g2SHpW0VNIESec1fpak8yWdKWlS6tk+kOI7oqHNRyX9OH3HVZLGDfS7NjiYYmukXSR9RdIXJXmnjhHCidg62S+A3STtt4M284Bf9R1ImgXMB46KiIP57V5k44DHIuItwBca2h8CfAI4FHgHcJqkVttWjQM2pPc3b9m0P3B5RBwAbAI+DXyZtDOwpJ2A/wp8LSI2RsTbI+KtwOLUFkkHUKybe3j6jm3AR3b0uzY5CHgauB24MyL+Iryi14jh0oR1rIgISZ8Cvpk6d+Mar0vajSKJXk6xASTAUcDNEfHL9Bn/ns7vStHDbvZO4JaI+HX6zGXAEbz6P/N/B/h3Wvt/UewVB/A14IyI+BtJz6SkvhfF1vbPpO94C/ANYE9gTnrf0cAhwH3pd92VIrF+fwe/a9992Jli3d2vA5+KiB/2E6d1KCdi62gR8R2grwywJ7Cy4fKZwBLgxTY+am/giSGEsi+wob8w+zm+GvhD4HcpesjFxYgHgd+T9GHgFKBvw8/rIuLPGz9I0rT0445+1wOA+4DXU/SkbYRxacJGqj2AE2hIcMl3gZMkvQFA0uvT+ZMoEl6zu4ETUl13N+D307lmJ5H+ILQwVdJh6edTgH9OP98CvJ9ix+DbUzy7S+rr2b/Ab3u3dwFzJU3ui1vSfxzgd+1zMHAPRfnjWkl79dPOOpR7xDZSvRH4XERsbXwmFRGPSDofWCFpG8UDrKeB3ShqstuJiPslfQX4cTp1dURsV5aQ9GmK+uyRkv6EokwxSdIHgdXAvwCfkfRlYA1wRfrsFyV9D3g2Ivp6qrOAJZL6dpb4k9R2jaT5wD+lmvJLwGeAp/r7XRscDPw4Iv5V0jnANyQdExEvtXkvLTPv0GE2AEkLgOURsbzh3HEUNd7lwHci4sAW79sJuB84KSIeG5ZgbURyj9hsYEspHpw1uh/Ypb83SJpJUcq4xUnYBuIesZlZZn5YZ2aWmROxmVlmTsRmZpk5EZuZZeZEbGaWmROxmVlmTsRmZpn9fwk/XsUVZOgiAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 360x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sub_hist(data['fr_number'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "#The approximate straightline form of the middle and tail parts (starting from k=15) of the histogram \n",
    "#evidences that the degree distribution approximately follows a power law. Let us calculate its exponent\n",
    "#using the formula (8.6) from the paragraph 8.4.1 of [1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "972609f98801483aa82754ae3a258bd3",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(IntProgress(value=0, max=616722), HTML(value='')))"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "2.4184901437551627\n",
      "0.0017524247495817025\n"
     ]
    }
   ],
   "source": [
    "k_min = 15\n",
    "\n",
    "N = data[data['fr_number'] >= k_min].shape[0]\n",
    "k_i_set = data[data['fr_number'] > k_min]['fr_number']\n",
    "\n",
    "list_to_sum = []\n",
    "\n",
    "for i in tqdm(k_i_set):\n",
    "    term = i / (k_min - 1/2)\n",
    "    list_to_sum.append(math.log(term))\n",
    "    \n",
    "Sum = sum(list_to_sum)\n",
    "answer = N / Sum + 1\n",
    "std = math.sqrt(N) / Sum\n",
    "\n",
    "print(answer)\n",
    "print(std)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#The exponent lies in the interval [2, 3], which is a common case for social graphs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Literature"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### [1] Newman M. Networks. – Oxford university press, 2018."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
